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ABSTRACT 

Many models of globular cluster formation assume the presence of cold dense clouds in early universe. Here 
we re-examine the Fall & Rees (1985) model for formation of proto-globular cluster clouds (PGCCs) via thermal 
instabilities in a protogalactic halo. We first argue, based on the previous study of two-dimensional numerical 
simulations of thermally unstable clouds in a stratified halo of galaxy clusters by Real et al. (1991), that under the 
protogalactic environments only nonlinear (8 ^> 1) density inhomogeneities can condense into PGCCs without being 
disrupted by the buoyancy-driven dynamical instabilities. We then carry out numerical simulations of the collapse 
of overdense clouds in one-dimensional spherical geometry, including self-gravity and radiative cooling down to T = 
10 4 K. Since imprinting of Jeans mass at 10 4 K is essential to this model, here we focus on the cases where external 
UV background radiation prevents the formation of H2 molecules and so prevent the cloud from cooling below 10 4 K. 
The quantitative results from these simulations can be summarized as follows: 1) Perturbations smaller than M m i n ~ 
(10 5 ' 6 MQ)(n? l /0.05 cm~ 3 )~ 2 cool isobarically, where rih is the unperturbed halo density, while perturbations larger 
than M max ~ (10 s MQ)(n^/0.05 cm" 3 )' 2 cool isochorically and thermal instabilities do not operate. On the other 
hand, intermediate size perturbations (M m - m < M pgcc < M max ) are compressed supersonically, accompanied by 
strong accretion shocks. 2) For supersonically collapsing clouds, the density compression factor after they cool 
to T c — 10 K range 10 2 ' 5 — 10 6 , while the isobaric compression factor is only 10 2 ' 5 . 3) Isobarically collapsed 
clouds (M < M min ) are too small to be gravitationally bound. For supersonically collapsing clouds, however, the 
Jeans mass can be reduced to as small as 10 5 5 Mq(tz7i/0.05 cm~ 3 ) -1 / 2 at the maximum compression owing to the 
increased density compression. 4) The density profile of simulated PGCCs can be approximated by a constant core 
with a halo of p oc r~ 2 rather than a singular isothermal sphere. 
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I. INTRODUCTION 

Globular clusters (GCs) are the oldest relics in our 
galaxy, presenting a cosmological challenge, as their 
ages often burst the allowed bounds of the expansion 
age. Since there are at least three distinct populations 
of GCs, the old halo and the middle-aged disk popula- 
tions of GCs in our Galaxy (Zinn 1985), and the young 
GCs in interacting galaxies (Ashman & Zepf 1992), we 
may well need more than one model to explain their 
formation and evolution. In the present study, we fo- 
cus on the formation of the halo GCs. Half of the to- 
tal mass of our galaxy appears to be contained within 
~ 50 kpc, whereas roughly half of the luminosity of the 
disk as well as half of the total number of halo GCs 
are distributed within the solar galactocentric distance 
(~ 8 kpc). So halo GCs clearly know about the dissipa- 
tion that followed galaxy formation and we must place 
their formation in the protogalactic context. This en- 
courages models that form halo clusters during the for- 
mation of the protogalaxy, such as thermal instabilities 
in the protogalactic halo (Fall & Rees 1985) or pre- 
galactic cloud collisions (Gunn 1980; Lake 1987; Kang 
et al. 1990; Kumai, Basu, & Fujimoto 1993). 



Many models of globular cluster formation assume 
the presence of proto-globular cluster clouds (PGCCs) 
in pressure equilibrium with hot halo gas (e.g., Gunn 
1980; Dopita & Smith 1986; Brown, Burkert & Tru- 
ran 1991; Kumai, Batsu, & Fujimoto 1993). The 
model by Fall & Rees (1985, FR85 hereafter) has 
been most widely adopted to explain the formation 
of such PGCCs: the thermal instability drives catas- 
trophic cooling of over-dense regions in the hot halos 
of protogalaxies and condenses them to dense clouds. 
In the present study, adopting FR85 model, we fol- 
low the thermal and dynamical evolution of over-dense 
clouds by numerical simulations in order to examine 
more quantitatively the physical properties of resulting 
PGCCs. 

Although the thermal instability was extensively in- 
vestigated in several areas of astrophysics, in recent 
years it has received much attention especially for X- 
ray cluster cooling flows: (cf. Balbus 1986; Malagoli, 
Rosner & Bodo 1987; Balbus & Soker 1989; Malagoli, 
Rosner, & Fryxell 1990). According to the results from 
numerical simulations in one-dimensional (ID) geome- 
tries (David, Bregman & Seab 1988 [DBS88 hereafter]; 
Brinkmann, Massaglia & Miiller 1990 [BMM90]), the 
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over-density undergoes a quasi-static compression in 
near pressure equilibrium when a cloud is small enough 
to adjust to pressure change faster than it cools, while 
larger clouds may undergo a "supersonic compres- 
sion stage" , leading to a density increase much higher 
than what is expected from the isobaric compression. 
The plane-parallel calculations by DBS88 showed that 
the density contrast increases to the "isobaric" ratio 
of the temperature of the hot background medium 
(Th ~ 10 7 — 10 8 K) and the minimum temperature 
(T c ~ 10 4 K) where the cooling becomes inefficient. 
On the other hand, the spherically symmetric calcula- 
tions of BMM90 showed that the density contrasts can 
be more than three orders of magnitude higher than 
those found in the plane-parallel calculations of DBS88. 
The difference between two geometries is mainly due to 
the fact that the gravitational acceleration increases as 
GM/R 2 when the spherical cloud collapses, while it 
remains the same (independent of R) when the plane 
sheet compresses. Also isotropic spherical compression 
(focusing) leads to much higher density accelerating 
the cooling and rise in density. However, 2D simula- 
tions of the collapse of an elongated blob by BMM90 
showed a shape instability to a pancake-shape leading 
to evolution best represented by the ID plane-parallel 
case. Thus geometry is clearly important in the non- 
linear evolution. On the other hand, 2D simulations by 
Vazquez-Scmadcni, Gazol, & Scalo (2000) showed that 
the condensation of fluctuating density field via ther- 
mal instability leads to a network of transient filaments 
which later accrete onto isolated round blobs. These 
blobs form at intersections between filaments just like 
galaxy clusters in large-scale structure formation in cos- 
mological simulations. This implies that the ID spheri- 
cal simulation can actually represent a generic collapse 
of the highest density peaks in the randomly fluctuat- 
ing density field in 3D. 

Numerical simulations in a 2D stratified background 
(Yoshida, Hattori, & Habe 1991; Reale et al. 1991 and 
references therein), on the other hand, have shown that 
the dynamic instabilities, such as Rayleigh- Taylor and 
Kevin-Helmholtz instabilities due to buoyancy- driven 
oscillations can disrupt the thermally unstable bubbles 
if the timescale for buoyancy oscillation (Brunt- Vasala, 
period) is shorter than the local radiative cooling time. 
Thus it is believed that the classic theory of thermal 
instability can not explain the effective mass depletion 
from the cooling flows of galactic clusters. We will ex- 
amine the role of such dynamic instabilities for the for- 
mation of PGCCs in the protogalactic halo environ- 
ment in the next section. The details of our models are 
given also in this section. The main results are in §3. 
Conclusions are given in §4. 

II. MODEL CALCULATIONS 



(a) Protogalactic Halo Model 

We adopt an idealized model for the protogalactic 
halo gas as in FR85; an isothermal sphere based on 
the Milky Way's properties, with a circular velocity 
V c = 220 km s _1 and corresponding gas temperature 
T h = 1.7 x 10 6 K. 

Many previous studies used a gas distribution spec- 
ified by the condition that the free fall time to the cen- 
ter of the protogalaxy equals the cooling time, yield- 
ing the halo gas density at a distance R g from the 
center p 1 ~ 1.7 x 10~ 25 g cnr 3 (i? 3 /10 kpc)" 1 . The 
use of pi has lead previous workers to consider halo 
densities of n h = 0.1 cm" 3 (FR85; Kang et al. 1990) 
or even rih — 1 cm~ 3 (Victri & Pesce 1995). For a 
primordial gas of H and He with an assumed ratio 
n(He)/n(H) = 0.1, the halo gas mass density is ph = 
(2.34 x 10~ 24 g cm~ 3 )n/j. For the isothermal sphere 
in hydrostatic equilibrium, however, the gas density is 
given by p 2 = 5.8 x 10~ 25 g cmT 3 f B {R g /lQ kpc)~ 2 , 
where fs is the ratio of the baryonic to total mat- 
ter density. Realistic values should be 0.02 <J rih <, 
0.07 cm -3 for R g ~ 10 kpc (assuming fs ~ 0.1) with 
some degree of central concentration. A value as high 
as 7j/j = 1 cm~ 3 is greater than the total density of 
baryonic and dark matter. Thus we take 0.05 cm~ 3 as 
a fiducial value and express the halo gas density as 



n h 



(0.05 cm" 3 )^" 9 



10 ! 



(1) 



where i?io is R g in units of 10 kpc and 1 <; q <; 2. 

Although the density irregularities in protogalaxies 
likely result from a superposed spectrum of perturba- 
tions, we consider, for simplicity, isolated spherically 
symmetric clouds resulted from sinusoidal perturba- 
tions. So the initial gas density of the cloud n c i ou d 
decreases gradually with radius r from its center to its 
edge R c as follows: 



TTcloud Tlh 



[..() + rfcosl 



for r < R c 



(2) 



where S is the amplitude of the initial density pertur- 
bation and R c is the cloud radius. 

The initial temperature throughout the cloud is set 
by the condition of hydrostatic equilibrium. To start in 
the nearly linear regime, we always begin with 5 ~ 1, 
unless otherwise noted. Initial parameters of the model 
clouds are summarized in Table 1. 

(b) Time Scales 

The free fall time tP^ is the time it takes a particle 
dropped from R g to reach the center of the protogalaxy, 
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Table 1. Initial Parameters for Model Clouds 



Model 


n h 


8 


$>c/ f>cool 


R c 


M c 






( cm d ) 






(PC) 




flow motion b 


SMI 


0.1 


1 


0.042 


4.20X10 1 


1.4xl0 3 


isobaric 


SM2 


0.1 


1 


0.092 


9.24X10 1 


1.5xl0 4 


isobaric 


SM3 


0.1 


1 


0.21 


2.10xl0 2 


1.8xl0 5 


isobaric 


SM4 


0.1 


1 


0.42 


4.20xl0 2 


1.4xl0 6 


supersonic 


SM5 


0.1 


1 


0.63 


6.30xl0 2 


4.9xl0 6 


supersonic 


SM6 


0.1 


1 


1.1 


1.05xl0 3 


2.3xl0 7 


supersonic 


SM7 


0.1 


1 


2.1 


2.10xl0 3 


1.8xl0 8 


isochoric 


SM8 


0.1 


5 


2.1 


2.10xl0 3 


3.2xl0 8 


isoclioric 


SM9 


0.1 


0.5 


0.42 


4.20xl0 2 


1.3xl0 6 


isochoric 


SM11 


1.0 


1 


0.17 


1.71X10 1 


9.8xl0 2 


isobaric 


SM12 


1.0 


1 


0.42 


4.22X10 1 


1.4xl0 4 


supersonic 


SM13 


1.0 


1 


0.68 


6.83X10 1 


6.2xl0 4 


supersonic 


SM14 


1.0 


1 


1.4 


1.37xl0 2 


4.9xl0 5 


supersonic 


SM15 


1.0 


1 


2.7 


2.73xl0 2 


4.0xl0 6 


isochoric 


SM16 


1.0 


1 


5.5 


5.46xl0 2 


3.2xl0 7 


isochoric 


SM21 


0.01 


1 


0.10 


1.02xl0 3 


2.1xlO b 


isobaric 


SM22 


0.01 


1 


0.12 


1.23xl0 3 


8.6xl0 6 


isobaric 


SM23 


0.01 


1 


0.21 


2.05xl0 3 


1.7xl0 7 


supersonic 


SM24 


0.01 


1 


0.33 


3.28xl0 3 


6.9xl0 7 


supersonic 


SM25 


0.01 


1 


0.41 


4.10xl0 3 


1.3xl0 8 


supersonic 


SM26 


0.01 


1 


0.61 


6.14xl0 3 


4.5xl0 8 


supersonic 


SM27 


0.01 


1 


0.82 


8.19xl0 3 


l.lxlO 9 


isochoric 




log(T) (K) 

Fig. 1. — The equilibrium cooling rate (f(T) = 
A/n 2 ) for a primordial gas composed of H and He is 
given in units of 10 _23 erg cm 3 s" 1 . Dotted lines show 
dln(A/n 2 )/dlnT = 1 and dln(A/n 2 )/dlnT = 2. 



The free fall time of a uniform gas cloud including only 
its self-gravity is 

(4) 

where n is the cloud density and any dark matter con- 
tribution is ignored. 

The cooling time for the halo gas of Th is given by 

t d ,h = 1 - 5n * kTh « (4.0 x 10 7 yrsK , (5) 

where the cooling rate of A/n^ = 5.5xl0 _24 erg cm 3 s _1 
at Th = 1.7 x 10 6 and the halo density given in Eq. 
(1) are used to calculate the latter relation. For an 
isobaric perturbation of the density contrast 6, the 
cooling timescale of the perturbed gas becomes t c i. c « 
t c l,h(l + <5)~ 2 , considering that the cooling rate is al- 
most constant for temperatures within factor of two 
of T ~ 10 6 K. We define the cooling distance as the 
distance over which the sound wave of the hot gas can 
travel within the cooling timescale of the perturbed gas, 

Icooi « c h t cl , c « (8.1 x 10 3 pc)i?? (l + S): 2 (6) 

where Ch — 198 km s _1 is the sound speed of the hot 
gas of T h . 
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The ratio of the cooling time to the self-gravity time 
for initial clouds is 



r/ = 6.2 x IQ- 2 R% 2 



l + S 



-3/2 



(7) 



For the range of density considered here, 0.01 < nn < 
1 cm~ 3 , cooling takes place too fast for self-gravity to 
be important initially. However, the cooling becomes 
inefficient once the gas cools to 10 4 K, and the gravi- 
tational time scale decreases as the density increases. 
The density in the condensed clouds increases by a fac- 
tor of Xisob = {T h /T c )((i c / n h ) ~ 350 for an isobaric 
case, where p c and ph are the molecular weights of the 
cold and hot gas, respectively. The density enhance- 
ment is even higher for supersonic collapses, resulting 
in the cloud pressure higher than the background pres- 
sure. The typical additional enhancement in supersonic 
collapses is 10 2 — 10 3 5 (see the next section). Thus, self- 
gravity becomes dynamically important in late stages. 
Also note that self-gravity would be more important 

— 1/2 

for lower density environments, since r\ oc n h 

The growth of thermal instabilities in a stratified 
medium is affected by the buoyancy-driven dynamical 
instabilities (see below), so the local buoyancy oscilla- 
tion period is also important. For our halo model it 
can be estimated by the local scale height H divided 
by the sound speed Ch, 



H 

Ch 



(5.0 x 10 7 yrs)i?i , 



(8) 



since H ~ R g for the power-law density distribution as- 
sumed here. For linear initial perturbation {e.g., 5 « 
1), the cooling time is about the same as the buoyancy 
oscillation period (t c i, c ~ tb) and the cooling distance 
is about the size of scale height of the halo (l coo i ~ R g ). 
This means the global stability needs to be considered 
(Balbus & Soker 1989) for such perturbations. 

(c) Effects of Dynamical Instabilities 

When a bubble-shaped perturbation condenses due 
to radiative cooling, it falls along the direction of grav- 
ity and so it is subject to the Rayleigh- Taylor and 
Kevin-Hclmholtz instabilities. Nonlinear evolution of 
such perturbations in a stratified halo of X-ray clus- 
ters was followed via numerical simulations (Hattori & 
Habe 1990; Yoshida, Hattori, & Habe 1991; Reale et 
al. 1991). The bubble is disrupted by these instabili- 
ties on a buoyancy oscillation period if it does not cool 
and condense significantly. Thus the cooling timescale 
must be shorter than the buoyancy oscillation period 
in order to prevent the disruption of the bubble. That 
depends on three parameters in order of importance: 
the ratio of the local buoyancy oscillation period to the 
local cooling time, the initial density contrast of the 
perturbation, and the ratio of the perturbation radius 
to the local scale-height (Reale et al. 1991). We note 
here the first and second parameters are not completely 



independent, since the cooling time of the perturbed 
gas is dependent on the initial density contrast. 

Following Reale et al. (1991), wc define for our pro- 
togalactic halo model, 



tb_ 

tcl,c 



1.25(1 + 5) 2 R 



10 



(9) 



using Eq. (6). According to their 2D numerical simula- 
tions, for £ « 1 the perturbations lead to fragmentation 
due to the dynamical instabilities, but the fragments 
can be thermally unstable. For £ ~> 5, the perturba- 
tion cools and condenses so rapidly that the dynamical 
instabilities do not have time to disrupt the bubble. 
These numerical simulations find that thermal insta- 
bilities cannot be an effective mechanism for mass de- 
pletion from the X-ray cluster cooling flows, since typ- 
ically £ ~ 0.01 there. According to Eq. (9), on the 
other hand, £ ~ 5 for 5 ~ 1 in a protogalactic halo, 
so the galactic halo provides more favorable conditions 
for thermal instabilities. 

Thus the initial density contrast of the perturba- 
tions is the next important physical parameter. Here, 
we will take the generic view that the protogalaxy was 
formed from the hierarchical clustering of smaller sub- 
structures and contains some dumpiness and turbu- 
lences. We suggest only nonlinear (i.e. S > 1) inho- 
mogeneities in the protogalactic halo can condense due 
to thermal instabilities without being disrupted by the 
buoyancy-driven dynamical instabilities. 

(d) Numerical Method 

The gas dynamical equations including self-gravity 
and cooling are: 



§ + P Vu = 0, 
at 

-n = —VP + g, 
at p 

de 1 , N A 

— = --V(Pu) + gu 

dt p K ' y p 



(10) 

(11) 
(12) 



where e = (7 — l)p/p + (l/2)w 2 is the total energy of 
the gas per unit mass and g is the gravitational acceler- 
ation. We use the Piecewise Parabolic Method (PPM) 
(Colella & Woodward 1984) in one-dimensional spher- 
ical geometry and treat gravity and cooling as source 
terms. The gravitational acceleration is softened near 
the center according to 



GM(< r) 



max 



(13) 



where G is the gravitational constant, M(< r) is the 
mass within the radius r, and s ~ 0.04P C is the soft- 
ening length. In most previous numerical simulations 
of the thermal instability (e.g., , the references in §1), 
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self-gravity was ignored because the gravitational time 
scale was longer than the cooling time scale and / or the 
emphasis was on the effects of the thermal instability 
alone. In our models, the gravitational time scale is 
longer than the cooling time scale initially, but not at 
later times. Thus, we must follow the transition of the 
inward flow from being driven by background pressure 
at early times to self-gravity in the last stages of col- 
lapse. 

The simulations start at t — with the clouds at rest 
in pressure equilibrium and cease when the background 
gas has cooled to I0 4 K (see below). The standard mir- 
ror condition for the reflecting boundary at the center is 
used, while the flow is assumed to be continuous across 
the outer boundary. 

(e) Radiative Cooling Rate 

Our radiative cooling rate f(T) applies to a primor- 
dial gas in ionization equilibrium (c/. Shapiro & Kang 
1987). (The cooling rate per unit volume A = n 2 /(T) 
where f(T) is a function of temperature only.) The 
success of FR85 model is hinged on imprinting of the 
gravitationally unstable mass scale of PGCCs at 10 4 K 
which is similar to the characteristic mass scale of GCs. 
This requires that formation of H 2 molecules should be 
delayed longer than a free fall time of PGCCs. Kang 
et al. (1990) suggested that the globular cluster might 
have formed only during the phase when the proto- 
galaxies are bright in UV background radiation. Hence 
we assume a model in which UV background radiation 
prevents formation of H 2 , so that the cooling rate is 
zero below 10 4 K and the minimum temperature is set 
at 10 4 K. 

The criterion for the thermal instability is that 
the over-density continues to cool faster than its sur- 
roundings, requiring d In f(T)/d In T < 1 to destabi- 
lize isochoric condensations and d\nf{T)/d\nT < 2 
to destabilize isobaric condensations (Balbus 1986). 
Fig. 1 shows f(T) and includes dotted lines where 
dlnf(T)/dhiT = 1 and dhi f(T)/d\nT = 2. 

We ignore the poorly understood process of thermal 
conduction. With the classical Spitzcr (1979) conduc- 
tivity, thermal conduction may dominate over the cool- 
ing. If we assume that the heat flux is saturated when 
the particle mean free path is comparable to the length 
scale of the temperature variation, the conductivity is 
one or two orders of magnitude smaller than the Spitzer 
value (see Gray & Kilkenny 1980). Brinkmann, Mas- 
saglia & Miiller (1990) showed that flow quantities dif- 
fer <^3% when they added the reduced thermal con- 
duction to their simulations of hot gas in clusters of 
galaxies. 

III. Simulation Results 



(a) Classification by Cloud Size 

We calculated spherical collapses for a wide range of 
halo densities (0.01 cm~ 3 < rih < 1 cm~ 3 ) and cloud 
sizes (0.04 <J R c /l CO oi <, 2.5). Initial parameters of the 
model clouds are summarized in Table 1, where M c 
is the mass of the cloud contained within the initial 
perturbation. 

Radiative cooling of the cloud leads to three differ- 
ent types of flow motions depending on its size rela- 
tive to the cooling distance, as recognized in previous 
studies (FR85, DBS88). According to our simulations, 
they can be summarized as follows: 1) (R c /l C ooi) ^ 0.2, 
small isobaric compression regime, 2) 0.2 <^ (R c /lcooi) Ss 
1 — 1.5, intermediate supersonic compression regime, 
and 3) {R c /l C ooi) ^ 1 - 1.5, large isochoric cooling 
regime. While in small clouds the collapse remains 
isobaric and quasi-static, in large clouds the central 
gas cools almost isochorically and the thermal insta- 
bility fails. Examples of these limiting cases are plot- 
ted in Fig. 2. Spherical collapses of clouds with 
Rc/lcooi ~ 0.04 and R c /l C ooi ~ 2 in the background 
halo with n h = 0.1 cm" 3 and T h = 1.7 x 10 6 K arc 
shown. Although the transitions between types are ob- 
viously continuous and gradual, we classified our mod- 
els into "isobaric" , "supersonic" , and "isochoric" types 
according to the characteristics of flow motions in the 
simulations. The last column of Table 1 shows the 
types. 

The intermediate size clouds belong to the most 
interesting regime. In order to see the evolution of 
an intermediate size cloud and the geometrical ef- 
fects, spherical collapse as well as plane-parallel col- 
lapse were calculated for one model: Th = 1.7 x 10 6 K, 
rih = 0.1 cm~ 3 , l coo i « 10 3 pc, and the scale height 
of the plane slab, H, and the radius of the cloud, R c , 
are H = R c = 410 pc ~ 0.4/ coo ;. Fig. 3 shows the 
evolution at t/t c i iC — 1.3,1.5,1.7,1.9, and 2.2, where 
t c i,c = 5.0 x 10 6 yrs. In both geometries, the gas 
cools isobarically until T ~ 2.5 x 10 5 K, below which 
the cooling becomes exponentially rapid and out of 
pressure equilibrium. As the pressure drops in the 
middle, a compression wave moves in. The infalling 
flow accelerates up to ~ 0.5c;,., moves to the center 
and then bounces back as an accretion shock. In the 
plane-parallel case, the accretion shock smoothes out 
the pressure gradient leading to near equilibrium with 
the background pressure. In the spherical case, on the 
other hand, the pressure of the shocked gas increases 
beyond the initial value owing to the central focusing. 
Thus the density enhancement factor can be as high as 
10 6 (n h ~ 10 4 cm" 3 ), which is - 10 3 5 higher than the 
isobaric compression ratio. Although there is a large 
positive pressure gradient between the cloud core and 
the background gas, there exists a halo of supersonic 
accretion flow instead of the out-flow, because the high 
density and high pressure core is gravitationally bound. 

Fig. 4 shows the density profiles of the clouds of 
intermediate size (R c ~ 0.42/ COO ;) with different halo 
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Fig. 2. — Evolution of small isobaric (model SMI: R c = 0.04i coo j, left panels) and large 
isochoric (model SM7: R c = 2.1l coo i, right panels) perturbations in one-dimensional spheri- 
cal geometry. Number density, pressure, radial velocity and temperature are plotted against 
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is time in units of 10 6 yrs. The small cloud is compressed isobarically with subsonic flow 
motions (no shock), while the large cloud cools with only little bit of compression. 
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densities, SM25 (n h = 0.01 cm" 3 ) and SM4 (n h = 
0.1 cm -3 ). This suggests that the density profiles of 
PGCCs can be represented by a core with a halo of 
p oc r~ 2 . Since the gravitational time scale decreases 
with the density as nT x l 2 while the cooling time scale 
decreases as n _1 , the relative importance of self-gravity 
increases in lower density models. For example, the 
model SM25 stops expanding after £ 6 ~ 215, while the 
cloud in SM4 expand outward slowly. But the expan- 
sion speed is about 0.7 km s _1 , much less than the 
sound speed of 10 4 K gas. 

Vietri & Pesce (1995) proposed that, instead of the 
quasi-static evolution into the two-phase medium, over- 
dense regions create an implosion shock that leads to 
H2 formation and cooling of the post-shock gas to 10 2 K 
(assuming no external UV radiation). However, we 
found that the implosion shock develops only if the ini- 
tial density perturbation is strongly nonlinear (S >> 1) 
or has a discontinuous top-hat distribution, both of 
which are unlikely in real density perturbations. Even 
in these cases, the implosion shock is nearly isothermal 
due to the rapid cooling of the post-shock gas. When 
the implosion shock rebounds at the center, it becomes 
an accretion shock. The post accretion shock gas still 
cools rapidly and condenses, so Vietri & Pesce's "small 
dispersing clouds" are an artifact of their methods. 

(b) Physical Properties of PGCCs 

According to our numerical calculations, only the 
intermediate size clouds in the supersonic compression 
regime, that is, 0.4 kpc <J i? c (n/ l /0.05 cm~ 3 ) <J 2 kpc 

and 10 5 - 6 M £ M c (n h /0M cm" 3 ) 2 £ 10 8 M , be- 
come PGCCs via thermal instability. Smaller clouds in 
the isobaric limit are not gravitationally bound (see in 
next subsection), while bigger clouds in the isochoric 
limit are too large to outrace the evolution of the back- 
ground gas to become distinct structures. The relevant 
cloud sizes increase for lower halo gas densities (so for 
greater R g ), since the cooling distance is inversely pro- 
portional to the gas density. But the formation time 
becomes longer at lower densities and the length scale 
increases to a scale on which other dynamical effects 
may suppress the density growth. 

In the transition from small to intermediate clouds, 
we see larger compressions owing to the supersonic in- 
fall with the compression factor increasing with the 
mass of the cloud (see Figs. 1 and 2). Thus the fi- 
nal radii of PGCCs, R P g CC , lie in the narrow range: 
Rpgcc {nh /0. 05 cm~ 3 ) ^30 — 60 pc. This owes to the 
dynamics of the supersonic infall whereby the shock in 
a larger cloud has a greater opportunity to accelerate 
and contributes more to its overall compression. This 
is extremely interesting, as the half mass radii of glob- 
ular clusters varies by only about a factor of ~ 3 over 
two decades in mass (c.f. Fall & Rees 1977). 

PGCCs form within 2t cLc « (2xl0 7 yrs) (0.05 cur 3 / 
while the background halo gas cools in a time scale 



of 2t cLh « (8 x 10 7 yrs)(0.05 cm" 3 /^). So PGCCs 
would lose the background pressure support in the same 
time scale, unless the halo gas is continuously heated 
by stellar winds and supernova explosion from Pop III 
stars or accretion shocks or shocks generated by merg- 
ing sub-clumps. Thus only self-binding intermediate 
clouds keep their identity, while clouds disperse. 

(c) Fragmentation of PGCCs 

Although our numerical simulations are limited to 
one dimension and cannot follow the thermal history 
of PGGCs below 10 4 K, we expect the dense core of 
PGCCs may fragment further owing to self-gravity. If 
the shape of PGCCs is more like a pancake, then the 
cloud would break up first by other dynamical instabili- 
ties into more or less spherical fragments and then only 
gravitationally unstable fragments undergo further col- 
lapses. Thus it is beyond our study to make the detail 
predictions on fragmentation processes. Here we at- 
tempt to make a rough estimate of the minimum mass 
scale of unstable PGCCs. 

For small isobaric clouds with M c < 10 5 6 Mq 
(n/j/0.05 cm~ 3 )~ 2 , the minimum unstable mass can be 
the "critical mass", M cr adopted by FR85, which was 
derived for an isothermal sphere confined by an exter- 
nal pressure (McCrea 1957), 

Mcr= i.i 8 ^y G -3/2 -1/2 (14) 

\pm H J 

For T c = 10 4 K and p h = 1.17 x 10- n dyne cm" 2 , 
M cr ~ 4.7 x 10 6 Mq > M sma u, so the small clouds are 
too small to become unstable gravitationally. 

For intermediate size clouds, as shown in Figs. 3-4 
(spherical cases), the density distribution of a typical 
PGCC can be approximated by a dense core with a 
halo of p oc r~ 2 . The flow in the core is almost static 
while the outer region is infalling supersonically. Thus 
the core of PGCCs can be approximated by an almost 
uniform static sphere. This is certainly different from 
the singular isothermal distribution (p oc r~ 2 all the 
way to r = 0), which is often assumed in the pre- 
vious studies mentioned in the introduction. In this 
case, PGCCs are confined by self-gravity and the back- 
ground pressure is negligible, so M cr given in Eq. (14) 
is not relevant. In a simple approach where the virial 
equilibrium is used (McCrea 1957) (e.g., 2K + CI w 0, 
where 2K = 3MkT/ (pmjj) is the thermal energy and 
O = — (3/5)GM 2 /R is the gravitational potential en- 
ergy for a uniform sphere) , the minimum unstable mass 
is given by 

M mr = 5.46 (j^^j ^ Pc- Xl \ (15) 

) ! where T c and p c are the temperature and density of 
condensed PGCCs. In fact this mass scale is almost the 
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same as Jeans mass of the isothermal uniform sphere 
(Spitzer 1979) defined as 

Thus we take this definition of Jeans mass as the 
minimum mass of unstable PGCCs as Mj ~ (5.9 x 

10 6 M )T 4 3/2 n~ 1/2 where T 4 = T c /10 4 K and n 2 = 
n c /(100 cm -3 ). Thus supersonic compression low- 
ers the Jeans mass by the square root of the ad- 
ditional compression factor. Owing to the spheri- 
cal convergence, the highest compression factor we 
have seen in our numerical simulation is 10 6 , up to 
~ 10 3 ' 5 times higher than the isobaric compression. 
For PGCCs considered here, 0.1 <J n 2 100, so the 
Jeans mass in the intermediate clouds is 10 5 5 Mq <, 
Mj(n h /0.05 cm -3 ) 1 / 2 <; 10 72 Mq. Considering that 
the critical star formation efficiency may range between 
0.1 and 0.5 in order for the resulting cluster to be self- 
bounded (Elmegreen et al. 1999), this Jean mass scale 
can be consistent with the observed GC mass distri- 
bution. If the evolution of more realistic clouds leads 
to irregular pancake-like shapes (§2.2, BMM90), how- 
ever, the density enhancement for real clouds would be 
smaller than what is found in the spherical simulations, 
but still larger than the isobaric compression. 

The use of shocks to lower the Jeans mass to a value 
that is closer to the characteristic mass of GCs was 
considered by Gunn (1980) in the context of chaotic 
protogalactic collapse and by Lake (1987) who invoked 
shocks created by cloud-cloud collisions. The cases con- 
sidered by Gunn and Lake both have so much shear 
that they would be unlikely to lead to bound structures, 
but shocks generic to the collapse of the intermediate 
mass clouds don't have this problem and still serve to 
reduce the Jeans mass to a value that is characteristic 
of GCs. 

IV. CONCLUSION 

We re-examined Fall & Rees (1985) model for the 
formation of PGCCs via thermal instabilities by numer- 
ical simulations of overdense clouds in a protogalactic 
halo. The key idea of this model is that characteristic 
mass scale of GCs, M c ~ 10 6 Mq, can be explained 
by imprinting of critical (Jeans) mass of gas cloud at 
10 4 K which has cooled from a hot halo gas in pres- 
sure equilibrium via thermal instability. If there were 
significant cooling below 10 4 K due to H2 and metals, 
the cloud will continue to cool and the imprinting is 
not possible (Kang et al 1990). Then the model fails. 
So we considered the cases where the radiative cool- 
ing is ineffective below 10 4 K, because the formation 
of H 2 molecules is delayed due to UV radiation from a 
central AGN or diffuse radiation from halo gas. 

The main conclusions are summarized as follows: 



1. Unlike X-ray cluster cooling flows, a protogalac- 
tic halo provides much more favorable conditions 
for thermal instability to operate. So PGCCs 
can form from marginally nonlinear {i.e. 5 J> 1) 
density inhomogeneities in a protogalactic halo 
via thermal instability. Linear perturbations are 
likely to be disrupted by the buoyancy-driven in- 
stabilities (Balbus & Soker 1989; Reale et al. 1991) 

2. The size of thermally unstable perturbations is 
determined by the cooling distance, l coo u over 
which a sound wave can travel in a cooling time. 
According to our numerical simulations in one- 
dimensional spherical geometry, the clouds of 
0.2/ coo ; <; R c <; l coo i can condense to form 
PGCCs via thermal instability. Smaller clouds 
(R c <; 0.2l coo i) cool only isobarically and become 
gravitationally unbound clouds, while larger cloud 
(i? c J> lcooi) cool isochorically and so do not con- 
dense to form any distinct structures. 

3. The thermally unstable clouds, which have mass 
10 5 - 5 M Q £ MpgccK/0.05 cm -3 ) 2 £ 10 8 Mq, 
are compressed supersonically, and the density 
enhancement of 10 3 — 10 6 results. This is much 
higher than what is expected from an isobaric 
compression. Considering the larger compression 
factor, we estimate the Jeans mass for PGCCs as 
10 5 - 5 Mq <, Mj(n h /0M cm" 3 ) 1 / 2 £ 10 7 - 2 Mq. 
The density distribution of simulated PGGCs can 
be approximated by an isothermal distribution 
with a constant core. 

Although we focus on the formation of old halo GCs 
in the Milky Way in this work, we note globular clus- 
ters exist in some extremely different systems. The 
largest of the dwarf sphcroidals, Fornax, Sagittarius, 
NGC 147, NGC 185 and NGC 205 all possess GCs. The 
specific frequency of clusters in these dwarf systems is 
as large as in bright ellipticals and larger than in spirals 
(Hodge 1988). If gas were in hydrostatic equilibrium in 
these parent galaxies, the temperature would be just 
10 4 - 10 5 K. Thus a different model other than ther- 
mal instability model for the formation of PGCCs or 
the formation of GCs themselves is necessary for these 
systems. 

Many models for GC formation have been suggested 
so far, but none of them seem to be complete or can ex- 
plain all the observed properties of GCs. One of the as- 
pects of FR85 model that has not been touched upon in 
the present work is the effect of the turbulent motions 
in the halo. Turbulent energy equivalent to the virial 
temperature should be injected at the largest scale and 
cascade down to smaller scales. According to 2D sim- 
ulations by Vazquez-Semadeni, Gazol, & Scalo (2000) 
the turbulence cannot suppress the thermal instability, 
because the time scale for the thermal instability re- 
mains shorter than the turbulent crossing time scale 
(^cioud/fturb) at smaller scales than the virialization- 
scale. This is because the turbulent velocity decreases 
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faster than the scale, so the turbulent crossing time 
increases at smaller scales. On the other hand, the 
turbulence can provide some heating to the halo gas. 
Since the halo gas cools radiatively in about 10 8 yrs, 
additional heating is necessary to maintain the halo 
temperature at Th for much longer than the cooling 
time scale in order to explain the observed age spread 
of GCs, typically a few Gyr. 

As discussed in §2, it is necessary to consider more 
realistic multi-dimensional simulations of thermally un- 
stable clouds in a protogalactic halo environment, which 
is currently under study. The present work, based 
on one-dimensional considerations, however, provides 
useful guidance for the up-coming multi-dimensional 
wotks. 
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